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Abstract 


In this paper, optimization of operating conditions for electrochemical energy sources is attempted using a dynamic optimization approach. This 
rigorous approach is demonstrated using a nonlinear diffusion equation that governs the mass transport limitation from/to the planar electrode. The 
resulting differential-algebraic model is solved using piecewise constant control vector iteration method that uses vectorized discharge current. 
The results obtained using the optimal control profile is compared with constant current discharge method and also with another trial and error 
approach that uses linear current control. It is found that the optimal control method achieves 12% more state of discharge against constant current 
discharge method with better energy efficiency and battery use. The proposed approach can be extended to complex systems such as operating 
electrochemical energy sources in hybrid environments that require proper control of energy distribution among the hybrid components. 
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1. Introduction 


Modeling and simulation issues in electrochemical energy 
sources such as batteries, fuel cells and super capacitors are 
gaining more attention these days. This is mainly because of the 
developments in hybrid power systems that have evolved as one 
of the major alternatives to using fossil fuels. An important topic 
in such hybrid environments is the power distribution among 
the hybrid components. Typically in many hybrid systems, an 
important component will be the electrochemical energy source, 
e.g., the batteries. Hence the optimization of operating condi- 
tions of electrochemical energy sources is important to control 
and share the power delivered by the batteries in hybrid environ- 
ment. Many previous investigations have attempted to optimize 
electrochemical energy sources using various approaches. Shep- 
herd derived a voltage—current and state-of-charge relationship 
for predicting discharge and charge curves particularly for lead- 
acid batteries [1]. Other forms of this relationship are derived [2] 
to relax the assumptions behind the Shepherd model. Another 
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relationship used to predict remaining capacity of a battery was 
given by Peukert’s empirical equation [3]. All these forms of 
simple relationships are not applicable to battery models in 
hybrid environments due to the fact that the battery can be 
studied only at limiting conditions under which these equations 
work. The main conclusion of a recent critical review paper is 
that the empirical equations cannot be used to predict state of 
charge accurately unless the electrochemical energy source is 
discharged under galvanostatic conditions at constant tempera- 
ture [4]. 

Current or potential pulsing is another attempt to optimize 
the charge—discharge characteristics of electrochemical energy 
sources. Purushothaman et al. [5] discussed the application of 
battery charging process using various current waveforms. The 
model consists of a simple linear diffusion equation with an 
analytical solution to study the mass transport limitations under 
different pulse currents. However, there are some disadvantages 
of pulse charging, in particular a decrease in charge amp-hour 
and energy efficiencies and an increase in cell electrolyte tem- 
perature. The constant current charge method results in good 
energy efficiency with no significant heat effects. It shows that 
an optimal operating policy with respect to the operation of elec- 
trochemical energy source that offers an advantage over constant 
current charging/discharging is unknown at this point of time. 
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Finding this optimal operating policy can provide the best way 
of utilizing electrochemical energy sources in any environment. 
This requires a more sophisticated way of keeping track of the 
remaining active material in the battery. 

In this paper, the problem of finding the best operating condi- 
tions is solved offline using open-loop control policy by defining 
the model equations in the form of a dynamic optimization 
problem. To illustrate the approach a planar electrode system 
is considered. Here, the control variable is discharge current 
density and the state variables are active species concentration 
distribution at node points across a discretized planar electrode. 
The method used to solve the dynamic optimization problem 
is the piecewise constant control vector parameterization. The 
optimal operating conditions obtained are compared with con- 
stant current discharging and the converged results are also 
verified using another solution methodology which involves the 
assumption of trial function for the control variable. The authors 
[6-8] have already adopted this approach for evaluating opti- 
mal operating conditions for electrochemical reactors used in 
electro-organic synthesis process. The rigorous dynamic opti- 
mization approach [9,10] will later be extended to combine 
concentration and potential gradients within a porous electrode 
and henceforth for the rigorous model governing the entire 
battery. 


2. Planar electrodes 


A mathematical model for planar electrodes that governs 
the mass transport to and from the electrode can be written as 
follows: 


c 

a © 
subject to the initial condition 
at t=0, c=co (2) 
and the boundary conditions 
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where D is the diffusion coefficient of diffusing species from the 
electrolyte, c the electrolyte concentration and i(t) is the applied 
discharge current density. The diffusion coefficient can be taken 
as constant in some cases, such as diffusion in dilute electrolyte 
solution. In other cases, such as diffusion in polymers, the dif- 
fusion coefficient depends very markedly on the concentration 
of diffusing substance, e.g., lithium-polymer. In this paper, the 
dependency of the diffusion coefficient is taken as 


D = Do (1 — 0.10) (5) 


The above model is based on the assumptions that: (1) the 
electrode is planar, (2) the electrode is completely charged at 
time t=0, i.e., 100% state-of-charge, (3) electrochemical reac- 
tion is taking place only at the surface of the electrode, (4) sym- 


metry at x=0 and (5) only the diffusion phenomenon is control- 
ling the process compared to electrochemical kinetics and other 
competing phenomena. The following dimensionless variables 
are defined to make the above model equations dimensionless: 
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Using Eq. (6), the dimensionless form of the model equation is 
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The boundary and initial conditions are 

at t=0, C=C (8) 
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The nonlinear model equations given by Eq. (7—10) can be 
solved by applying any discretization method in the spatial 
coordinate. Using finite differences in X with N number of inter- 
nal node points, the resulting differential algebraic form of the 
model equations are represented below: 
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where AX=1/(N+ 1) is the node spacing in X-direction. The 
initial conditions for differential equations (in Eq. (12)) are 
C;=1 at t=0 fori=1,...,N. The resulting differential-algebraic 
equations (DAE) represented by Eqs. (11)—(13) can be solved 
for any given discharge current density 6(t). The values of the 
parameters used for the simulation are length of the electrode 
L=8.2 x 107? m, diffusion coefficient Dp = 8.5 x 1078 m? s~!, 
initial electrolyte concentration co = 10 mol m7, constant dis- 
charge current i(t) = 10 A m~, equilibrium potential Ey = 4.1 V, 
discharge time tf = 530 s, number of electrons involved n = 1 and 
Faraday constant F = 96487 °C mol7!. 

The above set of model equations and set of parameters 
are simulated for the conventional case of constant current 
discharge, i.e., 6(t)=1. The resulting concentration variation 
across the planar electrode shows that the average concentration 
remains at a significant level of 62.11% state-of-discharge even 
after the end of discharge process. Obviously, the energy remain- 
ing inside the battery is critical and it should be utilized using an 
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effective alternate method of discharge other than an empirical 
relationship or pulse current. It is now clear that we need a profile 
based on the model equations and set of parameters, which can 
achieve a greater depth of discharge for the given time period. 


3. Linear current control 


Here, instead of using a constant current for discharge a 
trial function is assumed for 6(t). This is done by assum- 
ing that the discharge current varies linearly as a function of 
time. The parameters for any assumed policy can be evaluated 
using parameter estimation algorithm or using the simplified 
model equations in Ref. [11], developed by the authors. Using 
this approach, the linear current control is evaluated and bat- 
tery performance is investigated based on the model equations. 
It is found that a simple linear control can produce a depth- 
of-discharge equal to 69.11% which is 7% more than the 
depth-of-discharge achieved by the constant current control. 
Thus, this approach is interesting and further ensures existence 
of an optimal discharge current policy that needs to be evaluated 
using dynamic optimization technique. 


4. Optimal control policy 


The general statement for optimal control problem is as 
follows: 


min J[x(¢f)] (14) 
subject to the differential constraints 
x = f(x, u) (5) 


where x is a vector of state variables and u is a vector of control 
variables. J is the objective functional that is to be optimized 
and ¢ is the total process time. The specific objective function 
for planar electrode may be defined as maximize: 


Jo[5(t)] = 1 — Caye(te) (16) 


subject to the constraints represented as a system of DAE, Eqs. 
(11)-(13). Here, Jo is a function of t representing the perfor- 
mance index for dynamic optimization problem, Caye is the 
arithmetic average of concentrations obtained at all the node 
points (both internal and external nodes) and T+ is the duration 
of the discharge process. 

There are several optimization techniques [9,10] for the solu- 
tion of the above stated dynamic optimization problem, i.e., for 
solving Eq. (16) subject to the constraints Eqs. (11)—(13). Here, 
the method used for dynamic optimization of planar electrode is 
the piecewise constant control vector parameterization method. 
Since the problem considered does not have the nature of highly 
oscillatory control policy, a piecewise constant type of parame- 
terization is sufficient with the advantage of simple and efficient 
computation scheme. In this method the total time of discharge 
tr is discretized and the constraints are evaluated at each time 
step as follows: 

tf 
L=-— (17) 
p 


Initial conditions 
(concentrations, model 
parameters, etc.) 


Vectorized discharge 
current density, 5), 52, ... ,6p 


Solve the differential- 
algebraic equations 
(Eqs.[11]-[13]) using DAE 
solver (BESIRK) 


Evaluate optimization 
index, surface and average 
concentration at the end of 

operation 


Optimization algorithm 


Fig. 1. Computation scheme for evaluating optimal operating conditions. 


The computation starts at to where the initial conditions are 
known and ends at t, where p is the number of time steps. In 
each time step the control values to be used are determined by 
guesses that fall within the specified range for control variables. 
The strategy behind the random search in control domain is 
based on prior knowledge of the system. This search strategy 
can also enable quick convergence. Finally, the problem is also 
solved for different time steps, i.e., for p=10, 15, 20, etc. to 
ensure convergence in control domain. 

Fig. 1 shows the algorithm for developing the computer pro- 
gram. The model equations (Eq. 11—13) are solved using the 
DAE solvers DASSL in FORTRAN environment and BESIRK 
[12,13] in Maple® environment (for stability verification two 
different solvers are used). The number of internal node points 
used for the solution is 100 for better accuracy. The tolerance 
limits for the solver are also set high at 1078 for absolute tol- 
erance and 1075 for relative tolerance. The total number of 
equations needed to be solved for constant and linear current 
controls are 100 and 2 for the boundaries. But for the case of 
evaluating the optimal control profile, for each iteration, solu- 
tion of 102 differential-algebraic equations are needed. The final 
converged optimal profiles are obtained by using 10 runs with 
40 steps in time. This means that the solver has to stop at each 
time step and recalculate the state variables 40 times. Thus, 
one run requires the solution of 102 x 40 equations and for the 
optimal profile it is 102 x 400 equations. An executable code 
that calculates the electrode performance using different operat- 
ing conditions is available upon request from the corresponding 
author (VS). 


5. Discussion 


Fig. 2 shows the different operating profiles obtained using 
various approaches, i.e., constant, linear and optimal current 
controls. It also shows the corresponding potential variation 
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Fig. 2. Optimal discharge profile (a) and its corresponding potential profile (b) for a planar electrode. (i) Horizontal line, constant current; (ii) dotted line, linear 


current and (iii) continuous line with steps, optimal control. 


inside the planar electrode as a result of different controls. 
The electrolyte potential plot is made by using the Nernst 
relationship for surface concentration and electrolyte potential 
with an equilibrium potential assumed as 4.1 V. From these 
plots it is clear that a maximum load can be applied during 
the initial stages of discharge and then the capacity decreases 
gradually. The maximum level and the dynamic profile of 
gradual decrease in load sharing are obtained by keeping track 
of the active material remaining inside the battery. Absence 
of such an implied constraint prevents the constant current 
(or the other profiles) to achieve a greater depth of discharge 
in the same amount of time. Fig. 3 shows the concentration 
distribution obtained across the planar electrode at different 
time steps during discharge. Fig. 4 shows the average and sur- 
face concentration variations by following different operating 
modes in comparison with constant current control. It can be 
seen that the conventional constant current discharge method 
cannot utilize the active electrode material at the center of the 
electrode (x=0) due to transport limitations. But, a simple 
linear current control can yield better utilization than constant 
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control and thus reiterates the existence of an optimal discharge 
profile that can exploit the battery effectively and completely. 
Optimal control profile (curve c in Fig. 2) yields the best state 
of discharge of 74.93%. It is possible due to the dynamic current 
applied over the period of discharge that minimizes transport as 
well as kinetic limitations by having proper control over the 
amount of coulombs removed at any particular time step. More- 
over, the optimal control policy ensures the depletion of active 
material perpendicular to the electrode till the other end at x=0 
and discharges the electrode more efficiently. Figs. 3c and 4b 
show the concentration distribution for optimal operating condi- 
tions. From these results, it is evident that this method of keeping 
track of the availability of charged species cannot be properly 
achieved using heuristic controls like linear or constant cur- 
rent. The maximum states of discharge achieved using different 
methods are: constant current profile 62.11% and linear current 
69.11%. Table 1 shows the state of discharge obtained using 
different time steps in optimal control problem. It shows the con- 
vergence of optimization approach and its results in time domain. 
Fig. 5 shows the resulting power and energy distribution using all 
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Fig. 3. Concentration distribution inside the planar electrode at different times during discharge process: (a) constant current, (b) linear current and (c) optimal 


control. 
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Fig. 4. Average and surface concentration variations during discharge of a planar electrode using time-varying controls in comparison with constant current discharge. 
Continuous line denotes surface concentration and dotted line denotes average concentration: (a) linear current and (b) optimal control. 


Table 1 
State-of-discharge obtained using various time steps in control vector parame- 
terization approach 


Number of time steps, N Caye(te) SOD 

10 0.2953062825 70.46937% 
15 0.2712546656 72.87453% 
20 0.2653962833 73.46037% 
25 0.2530673100 74.69327% 
30 0.2507845663 74.92154% 
40 0.2507498113 74.92502% 


these controls. The maximum power and energy from a battery 
can safely be extracted by using optimal control profile followed 
by the linear and constant currents. Thus, the constant current 
discharge is not a suitable method for using electrochemical 
energy sources in varying load environments like hybrid power 
systems and some portable power applications with constant 
load fluctuations. For a given objective function, the optimal 
charge—discharge policy once evaluated using the physics based 


20 2.5 
2 
T 5 
3 ; 
g ae wo 
wo wn 
$ 10 2 
G A = 
Q Ko] 
N 7) 
= | = 
o] o 
E 5 £ 
[a] Q 


WN 
/ 


0 0.1 0.2 0.3 0.4 0.5 0.6 0. 
Dimensionless time, t 


Fig. 5. Power and energy distribution during discharge of planar electrode 
obtained using various controls. Continuous line, constant current; dotted line, 
linear current; continuous dark line, optimal control. 


model can provide a proper control over varying load hence 
improving performance and extending the life of the battery. 


6. Future directions 


This paper introduces a rigorous dynamic optimization algo- 
rithm for operating electrochemical energy sources particularly 
suitable in hybrid environments. To illustrate this approach a 
planar electrode system is considered and good performance 
improvement of up to 12% for state of discharge is demon- 
strated. The results are encouraging and the approach will be 
extended to porous electrodes governed by a rigorous battery 
model. Since temperature variation is one of the adverse effects 
in lithium-based batteries, it can be studied and optimized as 
another control variable using this approach. Such a study can 
provide new directions for the efficient operation of lithium- 
based batteries in hybrid cars. The nonlinear model considered 
for the demonstration is simple and it does not impose much 
computational burden. If we extend this approach to porous elec- 
trode or the entire cell model, more attention should also be paid 
for choosing an efficient solver in FORTRAN or other compu- 
tation environments. Also, the full battery model has 10 partial 
differential equations that cannot be directly used for dynamic 
optimization studies. This is because, if we choose 100 internal 
node points for discretization in both spatial and radial direc- 
tions then the model equations will be converted to a system 
of 4800 DAE [14]. During the process of optimization all these 
equations have to be solved at each time step and for each run 
covering the entire time domain. This may need the solution of 
at least 10 x 4800 x 100 DAE. 

To overcome this computation difficulty, an efficient simpli- 
fied battery model is being developed [11,14] to extend this 
optimization approach to a lithium-ion battery consisting of 
cathode/separator/anode. Incorporation of additional constraints 
related to power management is needed to set up an automated 
hybrid environment [15]. These and many other factors in imple- 
menting dynamic optimization procedure for battery models 


V. Boovaragavan, V.R. Subramanian / Journal of Power Sources 173 (2007) 1006-1011 1011 


will be identified and discussed in further publications using 
the approach introduced in this work. 
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